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Application of the magnetic field parallel to the plane of the graphene sheet leads to the formation 
of electron- and hole-like Fermi surfaces. Such situation is shown to be unstable with respect to 
the formation of an excitonic condensate even for an arbitrary weak magnetic field and interaction 
strength. At temperatures lower than the mean-field temperature the order parameter amplitude 
is formed. The order parameter itself is a (7(2) matrix allowing for the combined rotations in the 
spin and valley spaces. These rotations smoothly interpolate between site and bond centered spin 
density waves and spin flux states. The trigonal warping, short range interactions, and the three 
particle Umklapp processes freeze some degrees of freedom at temperatures much smaller than the 
mean-field transition temperature and make either Berezinskii-Kosterlitz-Thouless (driven either by 
vortices or half-vortices) or Ising type transitions possible. Strong logarithmic renormalization for 
the coupling constants of these terms by the Coulomb interaction are calculated within one-loop 
renormalization group. It is found that in the presence of the Coulomb interaction some short 
range interaction terms become much greater than one might expect from the naive dimensionality 
counting. 

PACS numbers: 73.63.-b,81.05.Uw, 72.15.Rn 



I. INTRODUCTION 

The fabrication of the graphene (graphite 
monolayers) and subsequent observation of the 
Integer Quantum HalP'^ in this layer produced a splash 
of theoretical and experimental activity. 

Though the transport properties of graphene arc con- 
trolled by the impurities, it is still worthwhile to un- 
derstand the phase diagram of the completely clean 
graphene. 

The main purpose of this paper is to point out that an 
application of magnetic field in graphene plane facilitates 
a spontaneous symmetry breaking. Though the dimcn- 
sionless Coulomb coupling /hv in graphene is large at 
large energies, it undergoes strong downward renormal- 
ization at small energies so that the analysis for weak 
magnetic fields can be safely carried out without resort- 
ing to any uncontrollable approximations^. 

To achieve this goal and understand the effects of the 
naively dimensionally irrelevant terms (such as trigonal 
warping, short range part of the interaction, Umklapp 
terms) we considered their logarithmic renormalization 
by the long-range Coulomb interaction and found some 
un-expectcd results. Earlier the effect of the Coulomb 
interaction was considered for the isotropic terms only^"^. 

The remainder of the paper is organized as follows. In 
Section II we discuss symmetries of a two-dimensional 
graphene sheet and write down its low energy Hamil- 
tonian. In Section III we describe physical reasons for 
the instability and suggest the order parameter. In Sec- 
tion IV we write down the Landau- Ginzburg free energy 
functional, discuss thermal fluctuations and describe the 



phase diagram. Section V is devoted to the analysis of 
the microscopies: we study the energy dependence of 
the effective coupling and renormalization of the lead- 
ing anisotropics in graphene. Section VI contains the 
summary and conclusions. Some auxiliary material is 
relegated into two Appendices. 

II. SYMMETRIES OF THE SYSTEM AND THE 
MODEL LOW ENERGY HAMILTONIAN. 



The purpose of this section is to write down the low en- 
ergy energy field theory to describe the electron-electron 
interaction in graphene. Our consideration will be based 
on the discrete symmetries of the lattice only and we will 
not appeal to any microscopic modcl^. 

The effective low energy field theory of graphene is 
constructed by the factorization of the original fcrmionic 
fields 5'o-(r;T),tT =t,i in terms of the oscillatory Bloch 
functions corresponding to the K, K' points of the Bril- 
louin zone, sec Fig. 1, 

«'<,(r;T) = V^^(r,T) *it(r). (2.1) 

where ij}c{v,T) is the four component fcrmionic field 
which can vary only over distances much larger than the 
lattice constant, and u{y) is the four-dimensional vector 
of the Bloch functions whose structure is described be- 
low. 

Consider the Bloch functions {w^ {r) k , {y) k) form- 
ing a basis for two-dimensional irreducible representa- 
tion of the wave- vector symmetry group C31,. (In the 



2 




4ti kr 



FIG. 1: a) The hexagonal lattice of graphene with the shortest 
translation vectors 01,2,3. Circles shows the positions of the 
maxima of the absolute values of Bloch functions u^'^(r). 
Their relative phases are shown, e = exp(i27r/3). b) The 
first Brillouin zone and the shortest translation vectors of the 
reciprocal lattice foi,2,3. Two non-equivalent Dirac cones are 
formed in the vicinity of points K and K' . 



tight-binding picture those Bloch functions are peaked on 
the corresponding sub-lattices, see Fig. 1). The overall 
point symmetry group is C^y and thus the wave-functions 
^K' — b^x] ' "H' = b^K] ^1^0 have to be included as 
points K' , K are connected to each other by C2 rotation 
and by time reversal symmetry. 

They are conveniently joined in a vector^° 



((^ 



A.„.B\ 
Ky^KjAB 



,B 



'k'/abIkk' 



(2.2) 



forming the basis of the four-dimensional irreducible rep- 
resentation of the planar symmetry group of graphene 
jo ^6^ii2i+j2a2 ! with the normalization condition 



dru{r) • u * (r) = 4|ai X 02]. 



(2.3) 



Hereinafter, f denotes integration within the unit cell. 

Thus, the fermionic field describing all relevant degrees 
of freedom, tp'^ = (V'f , V'i ) j has eight components. This 
eight-dimensional space is represented, see Eq. (2.2), as 
a direct product of the valley, {KK'), the sub-lattice, 
(AB), and the spin, (s), two-dimensional spaces. We 
will use standard 2x2 Pauli matrices fx,y,z, with the 
corresponding subscripts to parametrize all 8 x 8 matrices 
describing the Hamiltonian and the symmetry properties 

The partition function describing the low energy prop- 
erties of the interacting electrons in a clean graphene is 
given by {h= 1) 



l/T 



dr / dr 



OT 



(2.4) 

The Hamiltonian of the system has to satisfy all the 
discrete symmetries of the clean graphene, and to remain 
invariant with respect to transformations of the fields 
generated by the rotation C3, two reflections crj'^, and 



translations (i); 
C3 : V(r) - 

t : 



exp t/"- [ r X V + ^2.^ jj xp{r) (2.5a) 

i/>(x,2/)^S,A,-0(a;,-y) (2.5b) 
ip{x, y) t^kxxl}{-x, y) (2.5c) 



i/;(r) 



exp 



(2.5d) 



where 6**-^'* = ±27r/3. and we introduced the matrices 



(2.6) 



a = x,y, z. 

Continuous U{1) rotations in the spin space are given 

by 



U{1): t/>^exp i0,5,/2W> 



(2.7) 



where we choose z-direction of the spin to be along the 
magnetic field. 

Time reversal symmetry for the parametrization (2.2) 
acquires a natural form 

r : ^(r) r^^ ® ff ^' ® f^^^*(-r), B ^ -B (2.8) 

where B is the magnetic field acting in our case only on 
electron spin. 

Having listed the important symmetries of the prob- 
lem, we present the Hamiltonian in the form 

H = Hd + Hc + Hu. + Hsr + Hu . (2.9) 

The first term describes the Dirac-type kinetic energy 
and the Zeeman energy 



(re)V'^V-S?A-|-SV^S'^V- (2.10) 



= —iv 



where V = {dx,dy), S is defined in Eq. (2.6), the Bohr 
magneton and the 5-factor are included into the defini- 
tion of B, and Tc is the minimal linear scale present in 
the problem. As it was pointed out in Ref. 6 ( see also 
Refs. 7,8), the velocity v{rc) becomes scale dependent 
due to the Coulomb interaction 

H^^f. h,^ (^nr)^(r)) (V^t(ri)^(rr)) 
^ 2 J ^ |r-ri| 

whose strength cannot be renormalized as it is the 
only nonlocal term in the system. 

Though the terms described by Eqs. (2.10) and (2.11) 
are the most important ones on the dimensional grounds, 
they are not sufficient to define the problem completely 
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since their symmetries are much higher than of (2.5): 



Coo : ■4>{'r) ^ exp 



rx V + ^S, 



^Pir) (2.12a) 



CT„ : xP{x, y) ^ S.iAla;, -y) (2.12b) 
Uil) ® SU{2) ® 5f7(2) : -0(r) ^ J7 (a±, ^ ; m.a) V(i-), 
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exp 
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(2.12c) 

than it is allowed by Eqs. (2.5)-(2.8). Here 9' arc the 
continuous real variables and ni.2 are three-dimensional 
unit vectors. 

The terms lowering the symmetry of the Hamiltonian 
can appear both in the kinetic energy and in the interac- 
tion Hamiltonian. For instance, the trigonal warping of 
the one-electron spectrum is given by 



h.c. 



(2.13) 



where d+ = dx + idy, S+ = (S^; -|- iSj,)/2. Dimension- 
less coupling Xwifc) is of the order of unity at Tc of the 
order of the lattice constant and scales down at larger 
distances. 

In writing down the short-range interaction, one can 
neglect the effect of the Zeeman term on the scale of the 
order of the lattice constant. Thus the SU{2) invariance 
in the spin space must be preserved and ^ 



rcv{rc) 

+ E 

a—x,y,x 



E 

Q..(3—x.y,z 



Fo,p{r,){^^^^Kp4,f (2.14a) 



The symmetries (2.5) immediately yield the relation 



F — F —IF 

^ zx zy - 



F — F 

^ XX ^ yy 



Fx 



F,. 



J: 



2J 



J: 



yx — 
S,A. tE,A 



F^ 



J, 



-J 



S,A 



(2.14b) 



tS.A 



A reason for introducing F± , J± couplings will become 
clear later in Sec. V. As for the numerical values of the 
couplings, a very crude estimate at the scale Vc of the or- 
der of the lattice constants can be obtained by calculating 
the matrix elements of the bare interaction potential: 



Fal3 



\Av{rc) J J " J„c '"c|ai X a2||ri 
Pa/3(r) = (u(r)'l'SQA/3u(r)) , 



r2\ (2.15) 



and J^'^ = 0. (Obtaining finite values of J^'^ re- 
quires the virtual processes at least of the second or- 
der.) As Pa/3(r) contains only oscillatory components, 
see Eqs. (2.2) and (2.6), the integral in Eq. (2.15) is de- 
termined only by the distances of the order of the lattice 
constant. Thus, the parameters are extremely sensitive 
to the details at short distances and should be treated as 
entries for the low-energy theory. 

Though the warping and the short range interaction 
(2.13)-(2.14) lift most of the spurious symmetries (2.12), 
the extra continuous U{1) symmetry is still present cor- 
responding to Eq. (2.5d) with continuous 0*. It is related 
to conservation of quasimomentum which is violated only 
by the Umklapp processes. The lowest order term satis- 
fying the symmetries (2.5) has the form 



rlv{rc) 
6 



a/37 



(2.16) 

(^ts^A+^) (t/^tE^A+t/;) (t/>tE^A+t/>) + /i. 



where A_|- = [K^ + iAy)/2, Tap-y is symmetric with re- 
spect to permutations of the indices and 



^xyz — 0; Tzzz — 2T- + 

^ zxx ^ ziiit 'F— ^ F 4- 



(2.17) 



Finally, we notice that all the listed terms (2.10)-(2.16) 
remain invariant under the electron-hole transformation 



(2.18) 



This electron-hole correspondence will be important for 
the discussion of the instability arising in the in-plane 
magnetic field which we will discuss now^. 



III. PHYSICAL REASONS FOR THE 
INSTABILITY AND THE ORDER PARAMETER 

Having established the form of the Hamiltonian con- 
sistent with the symmetries of the lattice, we turn to 
the qualitative discussion of the instability and determine 
the target space of the order parameter. The symmetry 
arguments will allow us to do this without any actual 
calculation. 

Assume that no symmetries are broken at S = 0. 
Then, the low lying excitations are fermionic electron- 
and hole- like excitations with spin 1/2 and the disper- 
sion e(fc) = u(|fc|)|fc| as shown in Fig. 2 (a). The magnetic 
field parallel to the plane acts only on the spin and hence 
shifts the spectrum of the excitations making the creation 
of the Fermi seas for the electrons and holes energetically 



1 All the other short range spin rotational invariant interaction 2 rpj^^ leading irrelevant term lifting the e-h symmetry oc Vi/'^ 
terms can be reduced to those of Eqs. (2.14) by using the identity does not break any other interesting symmetries and will not be 

2(5o-io-2 5<T3 0-4 = <5(Ti (745^2 0-3 + Tctict4To-2CT3- important for our purposes. 
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favorable. Fig. 2 b). As the electrons and holes have the 
opposite charge, the Coulomb interaction makes them to 
attract each other. On the other hand, existence of the 
finite Fermi-surface leads to the Cooper-like instability 
first discovered in Ref. 11. This instability occurs even 
for an arbitrary weak interaction potential. 



b) 



hi fK.K' 




to the replacement A tj/S.tj'^. where C/ is given by 
Eq. (2.12c). Using Eqs. (3.1) and (2.12c), we obtain the 
most general form of the matrix order parameter 



A = Ao(r)t 



(3.2a) 



where Q is the 4x4 Hermitian matrix acting in spin and 
valley spaces subjected to the following constraints 



1 



- KK' 



1 



~AB 



(3.2b) 



The corresponding mean field single particle spectrum 
consists of four branches (see Fig. 2c): 



El^[e{k)±Bf + \^\ 



(3.2c) 



In terms of the angles in Eqs. (2.12c), the Q-matrix 
can be re- written as (a_, n2 — > a, n) 



Q = i^'^'(g) ( eif'' ) cos a 



-j<K' , '^2T" j sma 
(3. 2d) 

where ei^2 are two mutually orthogonal unit vectors in 
the plane perpendicular to the spin quantization axis 
ei = (cos 0s, sin 6*5, 0), e2 = (— sinf?^, cos6's, 0). 

Another way to parametrize Q from Eq. (3.2b) is to 
write 



Q 



V 
ft 



(3.2e) 



FIG. 2: The mechanism of instability in parallel magnetic 
field, B. a) The spectra of the one-particle excitation at S = 0; 
b) The shift of the spectra by finite B and the formation of 
the electron and hole Fermi seas; c) The attractive Coulomb 
interaction between electrons and holes leads to the instability 
towards a formation of the excitonic condensate^^ creating a 
gap in the one particle spectrum. 

As the result, gap Ao(T) is formed in the one-particle 
spectrum, Fig. 2 c). The resulting state is incompressible 
excitonic insulator with gapped charge excitations. The 
neutral excitations, however, are still quite interesting. 

As the electrons and holes can be paired with different 
phases and different valley indices can be involved, the 
resulting order parameter has a non-trivial matrix struc- 
ture which will be discussed now. As the exciton conden- 
sate is created by pairing of an electron and a hole with 
opposite momenta, the order parameter A has to be of 
the form (ipe'ipfi)- On the other hand the electron-hole 
transformation is defined in Eq. (2.18). Thus, we obtain 



A = (t/>®i/>t) = Ao(r)s,5, 



(3.1) 



Here A is 8 x 8 matrix acting in the space discussed after 
Eq. (2.3) and the Sz,^^; are defined in Eq. (2.6). 

If only symmetric terms (2.10) and (2.11) were present, 
the energy of the system would be invariant with respect 



where is a unitary 2x2 matrix in the valley space. 
Therefore, the order parameter is described by U{2) = 
SO{3) X ;7(1) matrices. 

Before writing down the effective action or the free en- 
ergy functional, it is better to explain the physical mean- 
ing of different angles in Eq. (3. 2d). To do so we will 
introduce the spin density and the "spin flux" , which in 
terms of the original (not smooth) fermionic fields have 
the form 



?(r) = i(*t(r)f^^,vl,,,(r)) 



(3.3) 



Jl J2 = l 



sm 



^(vl/t(r-fa,,)f^^,vl;„,(r + a,j), 



where the translation vectors ai_2,3 are shown on Fig. 1, 
and ji2 ^ ji - j2- 

Using Eqs. (2.1), (3.1), and (3.2), we find 



s(r) (X ezsina n^(|M^(r)p - |uB(r)p) 

+ 2?ixRcuA{r)uB{'r) + 2nj,Imu^(r)uB(r) 



(3.4a) 



The corresponding spin density configurations consistent 
with the phase factors of the Bloch function of Fig. 1 arc 
shown on Fig. 3 a)-c). Configuration of a = Tr/2, riz ~ 1 
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corresponds to the site centered spin density wave. It 
does not change the periodicity of the original lattice, so 
the Bragg peaks in the neutron scattering will remain at 
the same positions q = jibi + j2^2, and the ordering will 
affect only their internal structure. On the other hand, 
a — 7t/2, Tiz ~ 0, see Fig. 3 b.c). corresponds to the link 
centered density waves. It is easy to see that in that case 
the unit cell is tripled so that additional Bragg peaks at 
the positions q = ±K + jibi +^2^2 will emerge. Rotation 
in Ux-riy plane corresponds to a smooth transition be- 
tween Fig. 3 b,c) configurations. Such smooth rotation 
can be also understood as a continuous sliding of the su- 
perimposed spin-density wave with respect to the crystal 
lattice. 

The spin flux. Fig. 3 d) is maximal at a = 0, tt 



$(r) cx eiCosa(|MA(r)P - lusWp). 



(3.4b) 



Unfortunately, it is not coupled to any obvious physi- 
cal field which makes its direct experimental observation 
unlikely. 



IV. FREE ENERGY, THERMAL 
FLUCTUATIONS, AND PHASE DIAGRAM 

Order parameter (3.2) allows for smooth rotations be- 
tween all the states of Fig. 3, and thus is subjected to 
strong thermal fluctuations. Such fluctuations are gov- 
erned by the Landau free energy functional which we are 
about to describe. 

First we introduce a new definition of A, S-matrices 
Eq. (2.6), as 



An 



Sr, 



1 



KK' 



(4.1) 



where a = x,y, z, and A± = {Ax±iAy)/2. This modifica- 
tion is convenient because there are no gapless rotations 
possible in the sub-lattice space after the value of the 
mean-field order parameter is established. 

The only possible form for the free energy compatible 
with the symmetries of the system (2.5) is given by 



Z cx > ...m"" / VQcxp 



(4.2) 



where symbol . . . stands for a summation over 
the topological defects, sec subsection IV A, and fi de- 
notes a fugacity of such defects. A description of the 
thermodynamics of in terms of matrix Q subject to hard 
constraints (3.2b) is valid only at the distances larger 
than the spatial scale 



MF 



v{Rb) 

Ao(r)' 



(4.3) 



where the scale dependence of the velocity and the length 
Rb is defined in microscopic theory section, sec Sec. V. 




FIG. 3: a)-c) Sketclies of the spin density from Eq. (3.4a), for 
a — 7r/2, and a) Uz — 1 [N4 = 1 in notation of Eq. (4.5)], b) 
71a; = 1 {N2 — 1) and c) Tiy = 1 (A^3 = 1). Daslied arrowed 
lines show the primitive translation vectors of the resulting 
superstructure. Here, we choose 62 = (1,0,0), ei = (0,0,1), 
notice that the coordinate frame for the spin is rotated with 
respect to spatial frame (Axis x for the spin is the direction 
perpendicular to the plane and y is in plane axis perpendicular 
to the magnetic field.); d) Lines of the constant spin flux 
$(r) ■ ei (dashed lines) from Eq. (3.4b), a = (iVi = 1). 
Those lines can be also thought of as the spin current lines. 



The dominant term in the free energy. 



[Psiro) - PKiro)] 



-iTrS,Qd,Q 



(4.4a) 



[new notation is defined in Eq. (4.1)] has the symme- 
try which is higher than the symmetry of the original 
problem, see Eq. (2.12). Hereinafter, we will imply the 
summation over the repeating index v = x,y. 

To facilitate the further analysis, we rewrite Eq. (4.4a) 
using parametrisation (3.2e) 
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+ 2 [Psiro) - PKiro)] 



(4.4a') 



where is a unitary 2x2 matrix. Further investigation 
of free energy (4.4a) is postponed until subsection IV A. 

The remaining contributions, though may be small, are 
crucial because they remove artificial symmetries of the 
system. The following terms suppress SU (2) rotations in 
the valley space: 



= ^ T]a(ro)TrQkaQAa, rjx = Vy = V-L- 

07'q 

^ a—x.y,z 

(4.4b) 

]Fn = ^^ "^rQkpQkpQKqka. (4.4c) 



Finally, the term 
Ciro) 



4r§ 



Tr 



(4.4d) 



generated by the Umklapp processes, reduces the C/(l) 
rotations in the valley space, to discrete rotations (2.5d). 

A discussion of the role played by anisotropic terms 
(4.4b) - (4.4d) will be continued in subsection IV B. 



A. Isotropic part of the action and topological 
defects 

The sum over discrete topological defects can be re- 
placed by a path integral over continuous variables - dual 
fields. Introduction of dual fields become more transpar- 
ent when one uses the following paranictrisation of the 
order parameter (3.2e): 



y = e'^=5, g 



Ni H 

lN2 



lN2 

Ni- 



~N3 



(4.5) 



where g is an SU(2) matrix and J2i=i = 1- Then free 
energy density (4.4a') acquires the form 



i=l 



2 



(4.4a") 



So it appears that the U{\) and the SU{2) sectors of the 
theory are decoupled. This decoupling, however, breaks 
down when one takes into account the topological defects. 
As we shall demonstrate, the non-Abelian sector changes 
the selection rules for the vortices winding numbers. 

As usual for compact theories, the free energy density 
expression obtained in continuous limit has to be sup- 
plemented by topological defects to take into account 



a behavior of the order parameter in the vicinity of 
some singular points r^, where it vanishes: (3^(rj) = 
0; V{Yj)V'^[Yj) ~ 0. The absolute value of the order 
parameter is established at the distances of the order of 
^M_F, see Eq. (4.3), and at large distances, the defect is 
characterized by a contour integral around it. 

The topological defects are characterized by the wind- 
ing number 



1=^f dx,Tr{V+d,V) 



(4.6) 



The most standard approach would be to keep the g ma- 
trix in Eq. (4.5) single- valued and non-singular and create 
2tt- vortex configuration in C/(l) field 6*5: 



(x - Xj) ± i{y ~ yi) 



27r 



dx^dvOs 



(4.7) 



dx^Triy+d^V) = ±1, 



If only such excitations were present, U{1) and SU{2) 
sectors would remain decoupled and the standard 
Berczinskii-Kosterlitz-Thouless transition^^'^^ would oc- 
cur in the U{1) sector. 

There are, however, other configurations which are en- 
ergetically more profitable than those of (4.7). Indeed, 
consider the configuration of the form 



/2 



Jn-a/2 



)±Hy-yi 



^ {x-Xi)^ + {y-yip ' 

0; 



0' 



-in.-&/2 



(4.8) 



where n(r) is a smooth three-dimensional unit vector 

J2i=i = 1' ^^'^ Pauli matrices. 

Re- writing Eq. (4.8) in the form of Eq. (4.5), we obtain 
instead of Eq. (4.7) 



47r 



<j) dx^Tr{V+d,,V)±^. 



(4.9) 



i.e. configuration (4.8) corresponds to a tt or half-vortex 
in the U{1) sector glued with a half-vortex in the SU{2) 
-sector. Since the stiffness in the SU(2) sector vanishes 
at large distances, the change of sign of g in large defects 
does not require energy. So in the absence of anisotropy 
the main effect of the SU(2) sector is the change in se- 
lection rules of the vortices reflected in their topological 
invariant (4.6). 

To anticipate a role of the vortex and half-vortex con- 
figurations in the thermodynamic properties of the sys- 
tem we evaluate their energies. Substituting configura- 
tions (4.7) and (4.8) into Eqs. (4.4a') or (4.4a"), we find 



El = Trps In 



L 



IMF 



(4.10) 
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where L is the system size. As pK < Ps-, one can see 
that the half-vortices are in fact the main configurations 
responsible for the disordering of the U{\) sector, which 
we will in due course incorporate into the corresponding 
RG equations. 

To treat the contributions of the topological config- 
urations (4.7) and (4.9) systematically, we rewrite the 
partition function (4.2) (still neglecting the anisotropic 
parts) in the form 



Z (X / WVh cxp 



i / d\¥ 



PK 



[Ps - Pk] 
8T 



T 



T 



8ps 8p 



K 



Trd^h 



2pi 2/ii/2 sin7r|/i| 

— 2~ COS Zirhs 2 ~ cos TTAls ; 

''o '"o TT \h\ 



VV^ = 1; 



(4.11) 



where 2x2 matrix h is the dual field, t,^ y ^ are the Pauli 
matrices in the spin space, and /xi, /Xi/2, are fugacitics 
for the vortices and half- vortices respectively. Though it 
may appear that term Eq. (4.11) breaks the rotational 
symmetry, all the physical correlation functions calcu- 
lated with using Eq. (4.11) arc rotationally symmetric. 
Derivation of Eq. (4.11) is relegated to Appendix A 

Summing up the leading logarithmic series within the 
first loop renormalization group scheme, we find the fol- 
lowing equations for the corresponding fugacities: 



dp-i/2 
din ro 



T 



(PK + Ps) 



Ml/2- 



(4.12a) 



They are valid for M ^ 1 and also pk ^ T. 
The fugacity for the conventional vortices evolves as 



dpi 
din ro 



= 2- 



T 



Ml 



(4.12b) 



for /xi <C 1. Equations (4.12) are analogous for the simple 
estimate of the energy of the defects (4.10), however, they 
allow for a renormalization of the stiffness ps caused by 
the bound pairs of (half)vortices and (half)antivortices, 
and of the stiffness px caused by the bound pairs the 
half-vortices and anti-half-vortices as well as by thermal 
fluctuations of the order parameter. The latter renormal- 
ization, for pi/2,p ^ 1 and px ^ T, can be represented 
in the form 



dpK 
din To 

dps 
dlnro 



T 2 0(1/2)/ X 
P1/2PK {PK,Ps); 

TT 



(4.13a) 



p\,2fi'^'^\pK,Ps)-pWHps). (4.13b) 



Functions (3 in Eqs. (4.13) are difficult to obtain explic- 
itly for arbitrary stiffnesses. However, such forms will 
not be necessary for the further consideration. 



As the initial fugacities of the half- vortices, /ii/2 are 
not much smaller than the fugacities for the vortices /xi, 
the latter ones being less relevant can be neglected in 
the description of the phase transition. In the limit of 
Ml/2 1 one can also neglect the half- vortices in the 
renormalization of px- Then, from Eq. (4.13a) we obtain 
the following equations 



Pif(fo) = — In — ; 

TT ro 



S,K = iMF exp 



/ TrpiUip) 



(4.14) 

for (,MF ^ ^0 Ck- At ro ~ the pcrturbative renor- 
malization group is no longer valid, and the SU (2) sector 
enters into the strongly disordered regime. Hence length 
£,K lias a meaning of the correlation length of the order 
parameter. 

At ro 3> the non-Abclian stiffness vanishes: pK 
0. We can then use Eq. (4.12a) to obtain the value of the 
Kosterlitz jump at the phase transition: 



T^PsiTj 



KT 



-^.2*2- 



(4.15a) 



This value is modified in comparison with the pure U{\) 
model 



irpsiTxT - 0) 



(4.15b) 



due to the presence of half- vortices and the SU (2) sector 
being disordered. 

Assuming that Pi/2{(,mf) <C 1, we can also use 
Eq. (4.15a) to estimate Tkt'- 



t^Ps{Tkt) 



2T, 



KT 



2' + o[Mi/2(e 



MF 



(4.16) 



In Sec. V we will use Eq. (4.16) to relate Tkt to mean- 
field transition temperature Tmf- 

The main conclusions of these subsections arc: (i) spin 
U{1) and valley SU (2) rotations are coupled to each other 
in the presence of half-vortices, and (ii) the U{1) sector 
acquires the algebraic order at T < Tkt whereas the 
isotropic version of the valley rotations always remain 
disordered. 

The purpose of the next two subsections is to an- 
alyze how the disorder in the valley space is affected 
when the artificial symmetries are lifted by the leading 
anisotropics. 



B. Effect of weak anisotropies. 

Let us assume that the BKT transition in the U{1) sec- 
tor has already occurred so that the half-vortices are not 
relevant, and SU{2) sector is decoupled. (This assump- 
tion will be lifted in Sec. IV C.) Then the thermal fiuctu- 
ations lead to the renormalization of the anisotropics in 
the same fashion as the first term in the right-hand-side 
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of Eq. (4.13a). In the first loop approximation, wc find we have 



d hi ro 


2- 


2T 

TTPK_ 


da 


2 - 


8T ' 


din ro 


TIPK _ 


dC 


2- 


12T' 


din ro 


TTPK _ 



(4.17) 



Solving Eq. (4.17) with the help of Eq. (4.14), we obtain 



rg 

In^ 



^{r,) = n{iMF)[r^] \ —^] (4.18) 









i In^ 




( \ti^ 




Un/^ 



These equations are valid for ro ^ ^k. 

Thus, the effect of the anisotropics is determined by 
their value at distances of the order of ^k- If they are 
still small, i.e. 



C(6/f)<T 



S,MP 
^MP 
S.MP 



In 



K 



In 



In 



£,MP 
£,MF 



(4.19) 



MP 



then the SU{2) sector remains disordered and the con- 
clusions of the previous subsection remain qualitatively 
valid. If, however, one or several conditions (4.19) are 
violated, the anisotropies are important in determining 
the phase diagram of the system. 

The phase diagram obviously depends on the relative 
values of the coupling constants ij, n, Q and their signs. 
The latter are determined by the microscopic coupling 
constants (2.13) - (2.16), see also Sec. V. Since such 
couplings can not be established on a merely symmetry 
grounds, we will analyze the phase diagrams in the most 
possible general case and then construct the "physical" 
phase diagram in Sec. VE. 

Let us assume 



Vi^Mp) = Wai^Mp)] 

K,c 

Then, at the distance 
R* — ^MP 



1/2 > / £,MP 



K 



Imp 



(4.20) 



MP ) <;mp 



1/2 



(4.21) 



?7(i?,) ~ pk{R*), 



(4.22) 



so that the anisotropy becomes more important than 
the stiffness, and the perturbative treatment (4.17) is no 
longer valid. 

Instead, one has to identify the remaining soft-modes 
still compatible with the anisotropy potential (4.4b) and 
consider the fluctuations for those modes only. To 
achieve this goal, let us re-write Eq. (4.4b) using the 
parametrization (4.5): 



F 



1 

m 



(4.23) 



and 



and the isotropic part is given by Eq. (4.4a"), 

If there were no thermal fluctuations at all, the direc- 
tion of the order parameter would be obtained by mini- 
mizing the expression (4.23) with the results summarized 
on Fig. 4. 




FIG. 4: Phase diagram in the absence of the thermal fluctu- 
ation. For the pictorial representation of the states in terms 
of physical spins see Fig. 3. 

The regions (I) and (II) on Fig. 4, correspond to the 
Ising type anisotropy and the region (III) corresponds to 
the XY model. Hence deep inside the regions (/), (//) 
the system is ordered, whereas the region (///) is char- 
acterized by an algebraic order. At the lines separat- 
ing those regions the system possesses extra degeneracies 
which will be analyzed separately. 



1. Vicinity of the line r]± =0, rjz < 0. 

Near this line one can neglect the fluctuations of A^2,3 
and consider iVi_4 only. Then, the relevant part of the free 
energy, see Eqs. (4.4a") and (4.23), acquires the form 



PK_ 

2 



(4.24) 
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with iVi^4 constrained by N( + = 1. Representing 
A^i = COS0, 7V4 = sin0 and introducing dual field 8 to 
describe the vortices with a core size < i?*, wc find 



Zcx J V(j)Vecxp{-¥ii/T); 



T 2T^ 

. ^1 



2pK 



ida^edycj) (4.25) 



^2 J, cos 2(?i + ^ cos (27r6 



The derivation of Eq. (4.25) is completely analogous to 
the derivation of Eq. (4.f 1), see also Appendix. A. 

The first loop scaling equations for the anisotropy 77^ 
and the vortex fugacities /ii4 are 



In ro 
dr]± 



14 TrpK \ 
T I 
T 



= (2- 

2 - 



They have the obvious solutions 

'R 



(4.26) 



PiA {ro) = flu {R*) 



ro 



(4.27) 



If ri± (R*) = 0; the system undergoes Berezinskii- 
Kostcrlitz-Tliouless transition at Tbkt = t^PkI% from 
the high temperature vortex-dominated disordered state 
to a phase with power law correlations. 

At non-zero t]^ the ordered phase has a finite corre- 
lation length since below ■kpkI'^ the anisotropy r^x is a 
relevant perturbation. Therefore, phase </) is locked and 
the system is ordered with a finite correlation length ^14 
determined by the condition (C14) — T' so that 



14 



?7x(i?*) 



(4.28) 



From Eq. (4.28) one may conclude that the transition 
between phases (I) and (II) of Fig. 4 is a continuous one 
with smoothly varying correlation length exponent along 
the transition line. However, we will see shortly that this 
conclusion is an artifact of neglecting the higher order 
anisotropy (4.4c). 

At temperatures above 7rpif/2, both the anisotropy 
(order) and the vortices (disorder) are relevant. The cor- 
responding operators compete with each other since Q 
and (/) cannot be locked simultaneously. To estimate the 
location of the transition line and the length scale ^/ at 
which criticality becomes important, we require 



'7-L (6 



Mi4 (6) - 1, 



which yields, assuming as usual that /ii4(i?*) <C 1, 



P\^(R* 

ric{pK)^T[fiuiR*)]^'^^^ 



(4.29) 



At distances larger than ^j, the anisotropy is the 
largest scale in the problem, so that A'^i (iV4) becomes 
massive and the vortices provide the possibility for iVf 
{Ni) to change on the scale of the order of The re- 
sulting transition, therefore, involves only one soft mode 
and thus belongs to the Ising model universality class, see 
Appendix B, for more accurate calculation in the vicinity 
of TTpK = T. The corresponding correlation length ^ is, 
then, given by 



(4.30) 



r]JRJ 



<Nj^>=0 




K<0 

<Nj> = 
<N,>^0 



<Nj>4 
<N^> = 



K>0 



a) 



2T 




<N,>=0 



FIG. 5: Phase diagrams in the vicinity of the line ri± = 
0, rjz < 0. Lines A — b and A — c are always of the Ising 
type, a) Line A — d is (weak) first order phase transition line 
at K < 0. b) Lines A ~ d and A — d' axe second order phase 
transitions (Ising type) at «: > 0. At k = 0, the A — d line is 
a continuous phase transition with the varying indices given 
by Eq. (4.28). 

So far we ignored the higher order anisotropy term. In 
fact, those terms are important only in the vicinity of 
(I)-(II) transition line, where they change the order of 
the phase transition. 

To see this, we rewrite Eq. (4.4c) using parameteriza- 
tions (3.2e) and (4.5): 



7^- W ^ ^ 



-^Nyf^^cos^. (4.31) 
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Equation (4.31) should be added to Eq. (4.25). 

The couphng m Eq. (4.31) is renormalised by thermal 
fluctuations. The first loop RG equation and its solution 
is given by 



d In 7'o 
K(ro) = K(i?*) 



AT 

Ra, 



(4.32) 



At T:pK > 21" the quartic anisotropy k is a relevant per- 
turbation, so that the phase (j> becomes locked and the 
system is ordered in locally stable state even for ri± = 0, 
so that the correlation length ^14 is always limited by ^14, 

that is found from the condition k ( ^14 J ~ T. It yields 



a4<e 



14 



T 



(4.33) 



The length ^14 diverges when the approaching the multi- 
critical point npK = 2T, ri± = 0. At length scale large 
than ^14 the non-vanishing order parameter is found by 
minimisation of the expression 



F ^^±(114) 



T 



f2 rp 



cos 20 -I- sguK COS ^ 



(4.34) 



with respect to (j) which produces two second order (Ising 
type) phase transitions for k > 0, and the first order 
phase transition for k < [in the former case ^14 serves 
as pre-factor for the diverging correlation length similarly 
to Eq. (4.30)]. 

At 'KpK < 27" the quartic anisotropy becomes irrele- 
vant and we obtain from Eqs. (4.32) and (4.29) 



k{£,i) = k{R^) 



1 

/^14 



(4.35) 



i.e. it can not affect the position and the universality 
class of the Ising phase transitions. 

The resulting structure of the phase diagram in the 
vicinity of the line 77^ = 0, r/j < is summarised on 
Fig. 5. 



2. Vicinity of the lines rj^ > 0, 7?x = irjz . 

In the vicinity of the line 77^ > 0, 77^ = —77^ anisotropy 
(4.23) generates mass for the field N4 so that the free 
energy for the remaining soft modes, see Eqs. (4.4a") 
and (4.4d), is 



F : 



'2 — ^1 + Jl2 +iN3) ; 



(4.36) 



1; A7Ki?*) = 7y.(i?*) + r/i(i?0, 



where length scale i?» and the coupling constants are 
defined in Eqs. (4.21), (4.14), and (4.18). The vicinity 
of the line 77^ > 0, 77^ =77^ does not require a separate 
consideration as it is also described by the free energy 
density (4.36) after the replacement 



N^^Ni- 77x(i?,) ^ -r/i(ii'*). 



(4.37) 



At ro > R^, the scaling of the couplings in Eq. (4.36) 
is governed by the first loop rcnormalization group equa- 
tions 



dpK 


T 


din ro 


"2^ 




3T 


din Vq 


2 - 7, 

2ttpk_ 


dc 


21T ' 

2 


d In ro 


2npK _ 



A77; 



(4.38a) 
(4.38b) 
(4.38c) 



where we neglect the effect of the anisotropics on the 
rcnormalization of the stiffness which, however, will be 
sufficient for our purposes. 
The solution of Eq. (4.38a) is 



PKira) 



2-K ro ' 



Cff = R* cxp 



{2tt p{R^) 



(4.39a) 

which is consistent with the well known fact that the 
classical SU{2)/U{1) sigma model is always disordered 
by the thermal fluctuations with being the correlation 
length for these fluctuations. 

Anisotropics may lead to ordering. To estimate posi- 
tions of the phase transition lines, we solve Eqs. (4.38b) 
and (4.38c) with the help of Eq. (4.39a) and find 



A7;(ro) = A7;(i?,) 



ro_ 



C(ro) = C(i?*)(^ 



In 



ra \ 



(4.39b) 



In ^ 

^ R, , 



First, let us neglect the hexadic anisotropy, ^ = 0. 
Then, the anisotropy is important if |A77(ro ~ 2^^f)| > 
T. Equations (4.39b) and (4.39a) give us approximate 
positions of the phase transition lines 



A77x(i?,) 
T 

T 



2ttpk{R*) 
T 

2-n pk{R*) 



ai- 



T 



exp 



cxp 



T 

AnpKiR*) 



T 



(4.40) 



shown on Fig. 6. The numerical pre- factors ai^± of order 
of unity can not be obtained within a perturbative RG 
scheme. 

Line A77/ corresponds to the second order phase transi- 
tion from the disordered phase to the phase characterised 
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Ari 



<N,> = 



<N.> = 




FIG. 6: Phase diagrams in the vicinity of the line A?; = rjz + 
rj± =0, r]z > 0. Here A''+ = N2 + iN^. Lines Ari±^e are the 
Berezinskii-Kosterlitz-Thouless phase transitions. Line 77/ is 
the Ising type phase transition. Double line is the first order 
phase transition. Point C is the bi-critical point and point B 
is the critical end point. 



by (iVi) ^ 0. This phase transition belongs to the uni- 
versality class of the two dimensional Ising model. Addi- 
tional hexadical anisotropy, (finite but small Q can not 
affect this transition. 

Line Ar]± corresponds to the phase transition from the 
disordered phase to the phase characterised by algebraic 
correlations for iV2,3, i.e. to the Berezinskii-Kosterlitz- 
Thoulcss transition. Contrary to the Ising transition, the 
hexadic term may cause a further ordering and an addi- 
tional Bcrczinskii-Kostcrlitz-Thoulcss type transition be- 
tween the algebraically ordered to the long-range ordered 
phase. 

To describe those two transitions, we put A^i = 
in Eq. (4.36), introduce the field (f>, so that N2 = 
COS0, N3 — sin0, and the dual field 6 to describe the 
vortices. 

The singular part of the partition function is given by 



Zoc J P0I?0exp(-F23/T); 



T 2T ^ 2pK 



^ {da^ef + id^9dy(j) (4.41) 



c 



-^2^ cos cos {27:9) , 



where the coupling constant are described by Eqs. (4.14) 
and (4.39b) and R±^ < Rh found from the requirement 



Pk{R±) ^ A7]{R±_). This yields 

A?/_L(i?, 



R± R± 
^ In — = 



Ar^iR* 



1/2 



(4.42) 



Once again, the scaling of the vortex fugacity, (123, and 
of the hexadic anisotropy, C, can be determined from the 
first loop renormalization group equations 



dC 
din To 



2 - 
2- 



t^Pk 



T 
9T 



P23; 



(4.43a) 
(4.43b) 



In the fimit of /X23, C ^ 1; Eqs. (4.43) enable us to refine 
the definitions of a± (requiring '!TpKiR±) = 2T) and to 
find the boundary of the ordered phase r]e{R^) (requiring 
ttpk{R±) = 9T/2). With the help of Eqs. (4.42) and 
(4.39a), we find 



a± 



2.9; 



e 

2^73 ' 

/4e9\' 
V6-[7r^] 7?x«4.3*10V, 



(4.44) 



i.e. even though the boundaries of two phases have the 
same functional form, they are very well separated due 
to the numerical reasons. 

In the previous discussion, we imphed that C{R±) 1- 
This condition is clearly violated on the lines 775, rji 
if TipK /T 00, see the last of Eqs. (4.39b). Large power 
of the logarithm allows one to use the saddle point ex- 
pression and we obtain from condition C(-R±) ^ 1- 



7rp(i?, 
2T 



> 



21 



21 , T 
— In — — 

32 C(i?, 



1/2 



(4.45) 



At larger stiffness, the phase transition occurs by the 
locking of the order parameter due to the hexadic term 
and the fluctuations become unimportant. Minimising 
the potential part of free energy (4.36), we obtain a first 
order phase transition at r]{R^) ~ C(i?*). 

The only way of connecting the critical lines on the 
phase diagram, allowed by the symmetries of the system 
is shown in the inset of Fig. 6. 



C. Effect of strong anisotropies. 

In Sec. IV B we assumed that U{1) sector was ordered 
and the half-vortices were not important in ordering the 
SU{2) sector. It was justified by the assumption of the 
weak enough anisotropy so that the parameter ^, see 
Eq. (4.20), is constrained by f]{£_MF) <C T. In this case, 
the renormalization of the stiffness pn is so strong that 
Pk{R*) ^ Ps [see Eq. (4.21)] and separation of the scale 
was justified. We will see, however, in Sec. VD that for 
some ranges of the in-plane magnetic field the opposite 
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limits of the strong anisotropics, rj^^mp) > T arc more 
relevant. 

In this case the logarithmic renormalization of the stiff- 
ness pk, see Eq. (4.13a), and of the anisotropy parame- 
ters, see Eq. (4.17) are no longer strong and one has to 
investigate the soft modes in Eq. (4.23) already on the 
scale of the order of ^aif- Similarly to the approach of 
IV B only the lines of extra degeneracies requires further 
investigation. The difference is that U{1) sector and the 
half-vortex configurations have to be taken into account. 
To simplify further manipulations, we introduce the pa- 
rameter 



(4.46) 



All the results of the previous subsection correspond to 
the limiting case 7^0. 



1. The vicinity of line r]± =0, r/^ < 0. 

Let us generalise Eq. (4.25) by including the half- 
vortices. Substituting parametrisation (4.5) with A^2,3 = 
0, Ni — COS0, iV4 — sm<j> into Eq. (4.11), representing 
2h = diag [9 + hg; 9 — hg], we obtain 



Z cx y" V<i)D9V<P,V9s exp (-F^/T) : 



—TT COS [2nb 
Ri 



2pK 



RlT 



COS 20]- 



(4.47) 



'27T 

^cos (27r/is)] 



2pK 
Ml/2 

Ri 



COS (tt hs) COS {tt9) 



[to obtain fii4 term by generating vortices is completely 
analogous to Eq. (4.25)]. Here 7 is given by (4.46), /X14 
is a fugacity for creation of a vortex in (p field (we will 
call them "TVi — N4 vortices"), the vortices in 9s field [we 
will call them "vortices in the U{1) sector"] are governed 
by the fugacity ^1, and is the fugacity for the half- 
vortices. 

The first loop RG equations are analogous to 
Eqs. (4.12), and (4.26): 



dpi/2 
d In To 

dpiA 
din To 

dpi 
d In To 

dri± 
din ro 



2 - 



t^Pk /^l + 7 
AT 



7 



Pl/2] 



npK' 



r, ^PK \ 



(4.48) 



2- 



T 



We will see, that 7 = 1/4 is a special point where the 
A^i — A'4 vortices have the same scaling dimension as the 
half-vortices for the disordered valley sector. Another 
special point is 7 = 1/3 where the half- vortices become 
more relevant in the disordering the iVi — A'4 sector than 
the simple vortices. 

The phase diagram for 7 < 1/4 is shown on Fig. 7. The 
physics in the vicinity of the multi-critical point A is still 
determined by A^i — A'4 vortices whereas the half- vortices 
are irrelevant in this region. Inside the N1—N4 disordered 
region any deformation in (j) causes only a finite energy, 

{costt9) ^ 0. 

Therefore, the half- vortices in IJ(\) sector are allowed 
and the condition (4.15a) determines the value of the 
Kosterlitz jump. In fact, costt^ is proportional to the 
disorder operator p of the Ising model, see Ref. 23. We 
will use this analogy shortly. 

On the other hand, deep in the ordered sector, the half- 
vortices are confined and the disordering the t/(l) sector 
occurs due to the vortices and, thus, the Kosterlitz jump 
is determined by Eq. (4.15b). 

Only vicinities of the tricritical points I?i,2, 2 
the Fig. 7 a) deserve a special consideration. In their 
vicinity p\ term in (4.47) is irrelevant and p\4 and r\\_ 
terms conspire to produce an Ising critical point. There- 
fore at p\ii = we have two decoupled critical theories: 
the Ising one described by 0, 9 fields and the U(\) model 
described by {9s, hs) fields. At the Ising model critical 
point the operator cos 9 behaves as the disorder param- 
eter field of the Ising model and has scaling dimension 
1/8. Hence the corresponding RG equation for pif2 is 



dp 



1/2 



din To 



1 

2 - - 



T^Psira) 



AT 



Ml/2, 



(4.49) 



where we restored ps using Eq. (4.46). 

The scale invariance requirement for the fugacity pi/2 
gives the stiffness at the multi-critical points D1.2 



■np'sirQ 00) 
2Tc 



15 



< 4. 



(4.50) 



If Ps < Ps, the half- vortex fugacity pi/2 grows accord- 
ing to Eq. (4.49) as 



Pi/2{ro) _ f ro\ ^ 
Pi/2{R*) ^\R* 



(4.51) 



The grows (4.51) has to be stopped at the correlation 
length of the Ising transition 
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At the lengths larger than the Ising model becomes 
ordered, so the half-vortices are confined. However, 
their fusion produces usual vortices with the fugacity 
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algebraic 



C ' algebraic 
<N,>4 




FIG. 7: The phase diagram in the vicinity of the hne ryx = 
for different values of 7; a) 7 < 1/4; Zigzag line denotes either 
the first order phase transition or pair of the Ising transition, 
see Fig. 5. The asymptotics of the Berezinskii-Kosterlitz- 
Thouless lines are given by Eqs. (4.52) and (4.55) for points 
Di,2 and D'12 respectively, b) 1/4 < 7 < 4/15; The asymp- 
totics of the Berezinskii-Kosterlitz-Thouless lines are given by 
Eqs. (4.52) for points Di,2. c) 1/3 < 7 < 1/4. The difference 
between the slopes of the Ising and the Berezinskii-Kosterlitz- 
Thouless near multi-critical point B is given by Eqs. (4.52). 
Finite leads to the modification of the multi-critical point 
B. 



/ii ~ Mi/2(^j)- The requirement for the resulting pa- 
rameters of the remaining U{1) theory to be on the line 
of the Berezinskii-Kostcrlitz-Thouless transition leads to 



the estimate /i^^^(^/) — 1 or, see also Fig. 5 a), 
r]c Sp", 



In- 



Oil). 



For the spin stiffnesses, in the interval 

c c 2 ^Ps'^ 

Ps<Ps< Ps \ 



2tt 



(4.52) 



(4.53) 



the fugacity for the half-vortices (4.51) vanishes at 
00. If one tries to deviate from the line 77 = ijc towards the 
ordered side, those half- vortices fuse into vortices which 
are also irrelevant. 

Let us now consider the deviation from the tricriti- 
cal point towards the disordered side of the Ising sector. 
In this region cosO can be replaced by its average. As 
we have mentioned above, near the Ising model critical 
line cos 9 fi (the disorder operator of the Ising model). 
Hence we find 



cos7r6'cos7rft,. 



1^/^ costt/i. 



(4.54) 



The latter operator is irrelevant in the region (4.53) as 
well. Therefore, the line rj = rjc, < Ps < pf'^, is the 
transition line both for the valley and the spin sectors^^. 

At p > Pg'^, the half- vortices in the limit of the zero- 
fugacity pi/2 become irrelevant and U{1) sector is alge- 
braically ordered. The correction to vertical line can be 
evaluated using 



c,2 



Ml/2(C/) '^\V±-Vc\ 



1/8 



(4.55) 



This estimate is valid near points D[ 3. 

Case 7 > 1/4 needs further investigation as the critical 
points D'12 crosses the point A, see Fig. 7 b,c). For 
simplicity, we will neglect the quartic anisotropics: k ~ 0. 

Let us first consider 1/4 < 7 < 4/15. If r]j_{R^) = 0, 
the system undergoes the Berezinskii-Kosterlitz-Thouless 
transition and ttpk/^T = 1 where Ni — N4 vortices 
become relevant, (cos27r6') 7^ 0. At the same time 
(cos7r6') ^ and the remaining factor of the half- 
vortex operator becomes relevant, with the dimension- 
ality npK / i^'jT). As the result, the U{1) sector also 
becomes disordered. 

Initial steps in the consideration of ri± ^ are the 
same as in the derivation of Eqs. (4.29) and (4.30), and 
we obtain two lines of the Ising phase transitions. The 
scaling of the half-vortex fugacity on the Ising line is 
still governed by Eq. (4.49) so that the conclusions of 
Eqs. (4.50) and (4.52) remain intact, see Fig. 7 b). 

The peculiarity of 4/15 < 7 < 1/3 regime is that the 
half- vortex operator becomes relevant on the whole Ising 
line, whereas it is still irrelevant in the ordered region. As 
the result the points £'1^2 on Fig. 7 b) collapse, see Fig. 7 
c). The positions of the Berezinskii-Kosterlitz-Thouless 
lines still can be estimated using Eqs. (4.50) and (4.51). 

The most delicate case which we were not able to solve 
is 7 > 1/3. In this situation the half-vortices are al- 
ways more relevant than the vortices and the transitions 
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both in L'^(l) and in iVi — N4 sectors occur due to the 
half-vortices. The corresponding free energy obtained by 
keeping only half-vortices in Eq. (4.47) is 



+ ^ {d^&s)'" + ^ {d^Kf + td^h,dy0s (4.56) 



T7 



2jT 



2pK 



Ml/2 V± 

+ COS (tt hs) cos {nO) + cos 20. 

We do not know the critical property of this model of 
this model in the strong coupling lime. 



2. The vicinity of lines 7)z > 0, ri± = ±77^ . 

The analysis of those lines relies on the material of 
Sees. IV CI and IV B 2. Near the line 77^ > 0, A77 = 
|^±| —rjz'^ Vzj we have to generalise Eq. (4.36) [see also 
Eq. (4.37)] to include f/(l) sector. In this subsection we 
also neglect the hexadic anisotropy C = 0. The result- 
ing theory is obtained from Eq. (4.11) by putting extra 
constraint on the unitary matrix V 



= -V, Vi_ < 0; 
= ryVry, > 0. 



(4.57) 



One case case is mapped to the other by the substitution 
V iTyV and therefore only one of them has to be 
studied. 

The line A77/, see Eq. (4.40) is of the Ising type and 
all the analysis of Eqs. (4.49) - (4.55) is still valid. 

What remains is the vicinity of the Berezinskii- 
Kosterlitz-Thouless transition line Ari±, see Fig. 8. The 
theory describing phase transition in this case is, com- 
pare with Eq. (4.41), can be written in terms of dual 
fields only. 



Z cx y" VhsV0cxp (-Faa/T) : 



F23 T 2 



(4.58) 



M23(A??) 
i?2 



COS 2t:6 



Ml/2 (A??) 

i?2 



COSTt/i, COSTTf 



algeb. 



raic\ 




e'^^ algebraic 
<N,> = 



algebraic 




<Nj>4 

<^2,i> = 



FIG. 8: Phase diagrams in the vicinity of the line Arj = rj^ + 
^ ??z, riz > for finite value of 7, see Eq. (4.46), and in the 
absence of the hexadic anisotropics, C = 0. Here A''+ = A^2 -I- 



iA^3. Lines Arj 



X, 1,1/2 



are the Berezinskii-Kosterlitz-Thouless 



phase transitions. Line rji is the Ising type phase transition. 
The asymptotics of the lines Aryi 1/2 are given by Eqs. (4.52) 
and (4.55) for points D and D' respectively (those points 
are equivalent to the multi-critical points of Fig. 7 a) with 
the same notation. The transitions across the lines At^i^x 
are controlled by the unbinding of the vortices, and A771/2 is 
governed by the half-vortices. 



construct the phase diagrams in terms in the plane de- 
fined by anisotropies {riz{R*),ri±{R^,), see also Eqs. (4.20) 
and (4.21). Indeed, there can be no other singularities 
than those we have already considered, because for all 
the other regions of the phase diagrams the SU{2) sector 
is massive due to the anisotropies. Therefore, the char- 
acter of the singularities along the phase transition line 
remain the same. 

Bearing this in mind and expressing pk{R*) in terms of 
anisotropics using Eq. (4.22), we construct the "global" 
phase diagram shown on Fig. 9. The relation of this 
phase diagram to the physical coordinates B, T will be 
found in Sec. VE after the microscopic theory for the 
parameters of Landau free energy is built. 



where Ps ^ Pk{R**) because of the strong logarithmic 
renormalization at distances where the valley sector is al- 
most isotropic. It has a sequence of two phase transitions 
shown on Fig. 8 by lines A77X and A77i/2- 



D. Resulting phase diagram. 

In Sees. IV B, IV C, wc analysed of the vicinities of the 
degeneracies point and the vicinities of the degeneracies 
lines of Fig. 4. The results of this analysis enable us to 



V. LOGARITHMIC RENORMALIZATIONS 
AND MEAN FIELD TRANSITION 



The specifics of the problem in hand that it has three 
diferent intervel of the logarithmic rcormalizations: (i) 
energies larger than Zeeman splitting; (ii) energies larger 
than temperature but smaller than Zeeman splitting; (iii) 
classical renormalizations considered in the previous Sec- 
tion. As those renormalization are contributed by dif- 
ferent degrees of freedoms they have to be considered 
separately. The result will be the microscopic cxprcs- 
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Ising 

1 Berezinskii—Kosterlits—Thouless 
(vortices) 

1/2 Eerezinskii—Kosterlits—Thouless 
^^^=^^^= First order 

FIG. 9: The overall phase diagram for (a) decoupled valley 
and spin sectors, 7 — > 0; and (b) for 7 < 1/4. The positive 
quartic anisotropy, k > is assumed. The more detailed 
behaviour near the multi-critical points A,B,C,D are shown 
in more details on Figs. 5-8. The notation for those points 
here is consistent with that for Figs. 5-9. 

sions for the coupling constants in the free energy which 
will enable us to construct the physical phase diagram in 
Sec. VE. 



A. Logarithmic renormalizations at energies larger 
than Zeeman splitting. 

At such energies the Zeeman splitting is not impor- 
tant and can be considered perturbatively if necessary. 
We will also assume that the short range interaction and 
umklapp terms are not strong enough to lead to any re- 
construction in the state of the system at high energies, 
so they also can be considered within the perturbation 
theory, leading to a simple modification of the coupling 



constants which are not well-known anyway. 

The only terms which require the special attention are 
those related to the long-range Coulomb interaction. In- 
deed, a simple dimensional analysis of the Hamiltonian 
(2.10) - (2.11) points to logarithmic divergences in the 
simple perturbation theory (first identified in Ref. 6 for 
the three-dimensional gapless semiconductors). 

As usual, a summation of the leading logarithmic di- 
vergences is performed within the rcnormalization group 
scheme. However, the loop expansion would not be suit- 
able for the description of the graphcne as the dimension- 
less interaction strength e^/w(rc) is not small at distances 
Tc ~ a. 

Instead, the expansion in terms of iV - the number 
of independent fermion species entering the Hamiltonian 
(2.10) - (2.11) will be used. For the problem at hand, 
we have the valley and spin degeneracies so that N = 
4. The results, which will be obtained, indicate, that 
1/N corrections are quite small for = 4, so that 1/A^ 
expansion seems to be a reasonable approximation. 

To make the calculations compact, we will utilise the 
standard imaginary time diagrammatic technique^^. An- 
alytic expressions for the corresponding lines are given on 
Fig. 10. 

*^ =-G(e,k) 

e, k 

FIG. 10: Basic elements for the diagrammatic calculation for 
Hamiltonian (2.9)-(2.17), at distances < Rb = v{Rb)/B. 

Performing the rcnormalization group procedure, we 
change the smallest spatial scale in the problem from 
to the exponentially larger scale , It amounts to taking 
into account all the diagrams where the momenta q going 
through the Coulomb interaction propagator, belong to 
the region < \q\ < i/rf, This cut-off procedure, 

docs not violate the gauge invariance or other symmetries 
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of the problem. Then, we rescale fermionic fields as ■0 — s- 
(_1 + 5Z/2)il], ?/; ^ (1 + 6Z/2)-ip in order to keep the term 
^drip intact'^. 



a) 



d) 



e) 



sz = 



o 



A; = 



N 



6 O = 5Z X O 



16 + v{rc)^q'^ 



e = 




= 



/) 



□ 



5Z X 



□ 



FIG. 11: Diagrammatic representation for the leading in 1/N 
renormalization of the parameters of the high-symmetric part 
of the Hamiltonian, see Eq. (2.10). All the basic elements are 
defined on Fig. 10. The integration over the momentum q 
going through the wiggly line is restricted by 1/r^ < |g| < 
1/rf. 



and calculating the diagrams of Fig. 11 d), we obtain 
ding 



,1 — OAr/'j(5)l 

a lure Tr^iV 
where dimensionless function f-^ is given by 

arccos g 



(5.2) 



A.(5) = l-:r 



arccosh g 
9\f9 



,<?<!; 



(5.3) 



1 



, .9 > 1. 



It is easy to see that the function /„((?) is monotonously 
increasing and analytic for all g > 0, see Fig. 12. The 
asymptotic behaviour of this function is 



4 3 ^' 



1-^ 
2.9 



In 2.9 



oil 



.9 > 1- 



(5.3') 



Equations (5.2) - (5.3) were first obtained in Ref. 7 for 
= 1 where they are not applicable beyond the first 
term in the expansion (5.3') for g ^ 1. The validity 
of those formulas for TV ^ 1 was pointed out in Ref. 8, 
but the numerical coefficients here are different from the 
latter reference. 




FIG. 12: Plots of functions /u(<?) and fmig) entering into 
renormalization group equations. 



For the high-symmetry part of the Hamiltonian, see 
Eqs. (2.10) and (2.11), the only coupling which renor- 
malises is the velocity v{rc). This renormalization is a 
consequence of non-Lorentz-invariance of the Coulomb 
interaction. 

Introducing the dimensionless interaction strength 



9 (re) 



ne^N 
8w(re)' 



(5.1) 



^ In this scheme the scalar vertex and the Zeeman splitting term 
are not renormalised also, as a consequence of the gauge invari- 
ance. This can be checked by explicit calculation of the diagrams 
of Fig. 11 e,f). 



Solution of Eq. (5.2) with the help of the asymptotics 
(5.3') yields 



0.20 



(5.4) 



where TZ is the only relevant spatial scale generated by 
interaction. This scale can be estimated as 



In 



n 



In .9(a); .g(a) > 1 



nN 

May 



5(a) < 1- 



(5.5) 
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where a is the scale of the order of the lattice constant 
at which the continuous description becomes applicable. 

For the graphene sheets, the reported velocity, see e.g. 
Rcf. 17, is v{a) ~ 10® cm/s, we estimate g ~ 2^4 (uncer- 
tainty is associated with the dielectric properties of the 
substrate as well as uncertainty of linear scale at which 
the velocity is measured) and wc find from Eq. (5.5) 



7^~ 10-^ 



IQ-^a > a. 



(5.6) 



For all the further consideration we assume that the re- 
lation TZ ^ a is fulfilled. 

The Coulomb interaction strongly affects the scaling 
of the low-symmetry terms of the Hamiltonian. Let us 
start from the trigonal warping term (2.13). Calculating 
diagrams shown in Fig. 13a, we find 



d In A „ 



= -1 + 



7fwig); 



where negative monotonous analytic function 



(5.7) 



fwig) 



28 
15 

4 
9 



137r 



10 

^2 



6 



llTT 

' 2 arccos g 



~5 



X < 



2 arccosh g 
[ V9'-l ' 



5< 1; 
5> 1; 

(5.8) 



is also plotted on Fig. 12. 

The asymptotic behaviour of this function is 



bng 64g^ 



28 
15 



137r 



105 

10-8 In 2g 
? 



^ (5.8') 
Solution of Eq. (5.7) with the help of Eq. (5.8') yields 



Xwin) 

where 



1 



-5/8 



rc>n; 



(5.9) 



n 



1.57 



A»(a)(^)'" 



10"^ 10 



-5 



(5.10) 



Thus, we see that the Coulomb interaction tends to sup- 
press drastically the warping term making the energy sur- 
faces more and more isotropic^^. 

On the other hand, the Coulomb interaction leads to 
the enhancement on the short-range interaction terms in 



O =iZ X O 



<5 ! = 2dZ X > + 2 X ' + 2 X I ! + 2 X ^ 



+ 3 X I 




FIG. 13; Diagrammatic representation for the leading in 
renormalization of the parameters of the low-symmetric part 
of the Hamiltonian, see Eqs. (2.13), (2.14) and (2.16). All 
the basic elements are defined on Figs. 10, 11. The integra- 
tion over the momentum q going through the wiggly line is 
restricted by 1/r^ < \q\ < 1/rf. 



the Hamiltonian (2.14). Calculating diagrams of Fig. 13 
b) we find^ 

dlnFi'^ dlnjf-^ dlnJ^ 



dlnVc dlnrc 
dlnFl-^ d\nJ^ 



dhii 



d In d In 



f ( \ 



(5.11) 



where function fv{g) is defined in Eq. (5.3). 

Equations (5.11) can be easily solved with the help of 
Eq. (5.2) and we find 



f+'"(rc) _ J+'^(rc) _ J^jrc) _ a 
F^^^ia) J+'^(a) J- (a 
Fl'^ir,) J^(r,) afgia) 



F^' (a) J^{a) Tc \girc) 



(5.12a) 
(5.12b) 



Couplings in the Eq. (5.12a) are irrelevant. The inter- 
actions in Eq. (5.12b) arc strongly enhanced by the long- 
range Coulomb interaction. This enhancement becomes 
especially pronounced at intermediate distances Vc < 72.. 
Using Eq. (5.4) we find 



ff-^(rj _ J^(rc) 
Ff'^(a) 



^ Notice that the "mean-field analysis" of Ref. 5 of the excitonic 
instabilities at zero magnetic field corresponds to accounting of 
only the third diagram in the right-hand-side of Fig. 13 b) and 
thus is false even within 1/7V approximation. 
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i.e. naively dimensionally irrelevant couplings become 
weakly relevant. Though, such small value of indices is 
clearly beyond the accuracy of the 1 /N approximation for 
iV = 4, this formula indicates, however, that the effect 
of the short range interaction lowering the symmetry of 
the system is much stronger than it was thought before. 
Whether or not this enhancement may lead to instability 
at zero magnetic field requires further improvement of 
the renormalization group scheme, which is beyond the 
scope of the present paper. Here, we simply note that 
at large distances, g — > l/ln(rc), see Eq. (5.4), so that 
the coupling of Eq. (5.12b) becomes irrelevant again. It 
indicates that at zero magnetic field the excitonic insta- 
bility can occur only as a first order phase transition. In 
all subsequent consideration, we will assume that such 
transition does not occur. This assumption is in accord 
with all the experimental findings accumulated so far. 

For the further use, let us recast the answer (5.12) 
for the most important constants, in the form similar to 
(5.9). Because of estimate (5.6), the second digits in the 
exponents are not observable and will be omitted: 



<7^; 



Fl'^{TZ)^Fl'^{a) ~ 1. 



1]'; re>7^; 



(5.13) 



Finally, the Umklapp terms are also enhanced by the 
interaction. Calculating the contributions shown on 
Fig. 13 b), we obtain the renormalization group equa- 
tions 



dlnrc 
din Tr 



-3; 

-3 



64 



Mg). 



(5.14) 



Similarly to the (5.12), the solution of Eq. (5.14) is 



.9(a) 



T+{a) \rcj ' ^+(a) \rcj V5(^c), 

(5.15) 

In particular, at the intermediate distances Tc ^ TZ, we 
find 



3-- 



1.4 



i.e. the three particle Umklapp interaction remain ir- 
relevant though it is strongly enhanced by the Coulomb 
interaction. 

To conclude this subsection we notice that the short 
range interaction terms are vital in the consideration of 
the QHE ferromagnets^^'^^ and the effect of warping on 
the weak localization was considered in Ref. 20. The re- 
sults of this subsection indicates that the estimates done 
in those works arc hardly reliable. 



B. Logarithmic renormalizations at energies 
smaller than Zeeman spHtting: separation of the 
electron-hole and the Cooper channels. 

In the previous subsection we considered the Zeeman 
splitting as a perturbation. This is legitimate to do up to 
the spatial scale < Rb, where the length Rg is found 
from the equation 



Rb\ 

n 



BRt 



(5.16) 



The scale dependent velocity v is determined from 
Eqs. (5.1) and (5.2) and we highlighted that this depen- 
dence may include only one spatial scale TZ. To solve 
Eq. (5.16) and facilitate further discussion, we introduce 
the natural scale for the Zeeman splitting, So, according 
to 



-Bo 



v{rc ^ n) 
n 



(5.17) 



rough estimate for Bq is Bq ~ 10 ^ lO^K. Then the 
solution of Eq. (5.16) takes the universal form 



Rb = n.fB ( 



(5.18) 



where the function fsix) is the solution of the equation 

V Ub {x)] = xv{l) 
/s(a;)«l/a;i-25; x > 1; 

/M.)-^in(i)+i; x<l; (5.19) 

At > Zeeman splitting freezes some electronic 
degrees of freedom, and, on the other hand, gives rise to 
the finite density of state for the other ones, see Fig. 2 
b). 

To make use of such separation, we will include the 
Zeeman splitting in the denominator of the Green func- 
tion and decompose the result as 



G 



1 



G = V{n)- 



ie~v{RB)kJ:- BSz 
1 



ie - 

^ = B ~ v{RB)\k\; n - 



-SG; 
V{n) 



1 



k 



2B) 
(5.20) 



where 



V 



1 — n . Y,S, 



V; 



(5.21) 



is the projection operator to the branches of the elec- 
tron spectrum which may produce excitations with the 
energies smaller than B. 
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The 5G component of the Green function (5.20) does 
not have a resonant denominator and may be neglected. 
For tlic logarithmically divergent contributions only ^ < 
B are important so that the integration over the momen- 
tum can be replaced by 



{2nf 



dfi 



27r 2Trv^{RB) 



(5.22) 



Usual calculation of the second order correction to the 
interaction vertex shown on Fig. 14 a,b) reveals two loga- 
rithmically divergent contributions, which can be readily 
identified as electron- hole (a) and Cooper (b) channels. It 
is worthwhile to emphasise that those contributions are 
associated with the presence of the Fermi surface, and 
have the structure very different from the high-energy 
logarithmic terms of the previous subsection. 



= -g- Q = P{n) 



-£,-fc:T (i) 



be summed up in all orders of the perturbation theory. 
The fact, that the logarithmic divergence occurs only in 
one channel justifies the mean-field approximation which 
will be employed in the next subsection. 



C. Mean field transition. 

After the leading divergent series is identified, it can be 
summed within the standard mean field approximation^^, 
shown on Fig. 15a-c). Simple examination of the momen- 
tum q and frequencies lo transfered through the interac- 
tion wiggly line shows that all the momenta q < 
contribute almost equally into the formation of the or- 
der parameter, whereas the approximation w w is valid 
with the logarithmic accuracy. The finite order param- 
eter A changes the polarisation operator only for small 
q ~ A/v{Rb) ^ This modification may produce ef- 
fect only of the order of A^/S^ and it will be neglected, 
see Fig. 15 c-d).^6 

= -QA{e„,ky, e„ = 7rr(2n+l); 



d) 





111 



FIG. 14: Leading (a-d) and sub-leading (e) logarithmic di- 
vergences at Tc > Rb- Series (a,c) corresponds to the exci- 
tonic instability^^ in the electron-hole channel studied in the 
present series. The interaction in the Cooper channel (b,d) is 
repulsive and, therefore, renormalises to zero. 

To collect the leading logarithmic divergences, we look 
at the third order diagram, and find that only those 
corresponding to the ladder series, see Fig. 14 b,d) are 
proportional to the second power of the logarithm. The 
sub-leading terms, see Fig. 14 e), can be combined as a 
perturbative (in 1/A'^, or the interaction strength) renor- 
malization of the coefficients in the second order diagram 
and, therefore, can be neglected. 

Finally, it is easy to check that the signs of the dia- 
grams Fig. 14 a) and c) are the same, whereas, the signs 
of the diagrams Fig. 14 b) and d) are different from each 
other. The latter corresponds to the repulsive interac- 
tion in the Cooper channel, which, thus renormalises to 
zero. The former one describes the attractive interaction 
of electron and hole and leads to the excitonic instability. 

Therefore, the only relevant diagrammatic series is the 
ladder series in the electron-hole channel, which should 



6) A = -A(n) = 



o 



, Ao 1,2 \ N B 



d) 



LO, q 



NB 



8g{RL 



iT\q\RB + 8g(i?_B 



FIG. 15: Mean-field equations corresponding to the summa- 
tion of the most divergent series. Fig. 14a) , c) . 

The resulting mean field equations are obtained from 
Fig. 15 a-b) with the interaction propagator Fig. 15 c) 

using the approximation |pi — p| ~ 2/ Rb sin ^n; ni, valid 
once again with the logarithmic accuracy, and the inte- 
gration rule (5.22). It results in 



Bd(, ze„-?^, -HA(r?i) 



<L„=irT(2n+l 

X r{ni) 



B 2^v-^{Rb) el +e+ ^(ni)A(ni)2 

(5.23) 



where 



V{a) 



2tiv\Rb) 



89 



NB 

Vn 



^27rsin^ -t- : 



r da , , 
n = 7;— ^(a) cosna. 

J~7T ^TT 



(5.24) 
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Substituting A(n) of the form 

A(ri) = S^±, {Ao(T)^(ri) + Ai(T) [l - } 

(5.25) 

into Eq. (5.23), we find 

and the self-consistency relation involving Ao(r) only. 
At T = the gap is given by 



^ {Vo + Vi)B 2B 



(5.26) 



where the constants Vb,i are given by Eq. (5.24). 

Substituting Eq. (5.24) into Eq. (5.26), we obtain after 
the integration 



Ao(T = 0) = 2Bexp 



2N 



(5.27) 



ig(.RE 



where /a (2:) < 1 is a dimensionless function of the scale 
dependent interaction strength (5.1), and Rb is defined 
by Eq. (5.16). 

The explicit expression for this function is 



/a (a:) = < 



■ arccosh x 



1 TT 

X 2x^ 



x^ 1 
arccos x 



X 2x^ 



x>l 

X <1 
(5.28) 



and its asymptotic behaviour is given by 



4.T , o> 



a; < 1 



4 , 2a; 2 / 1 , 
— In — I' — + { — ], x:$>l 

TTX e x'' \ X'^ 



(5.28') 



To write down the explicit expression for the zero- 
temperature gap, we use the scale of the magnetic field 
introduced in Eq. (5.17). If the magnetic field is weak, 
S <C So, then Rb S> TZ, and the effective interaction is 
also weak, g{RB) ^ 1- Using asymptotics Eqs. (5.4) and 
(5.28') in Eq. (5.27), we find 



Ao(T = 0) w 2Scxp 



(5.29a) 



41n + 

For the strong magnetic field, B ^ Bq, the result reads 



Ao(r = 0)«26cxp<^ -2N 



3 \ B 



(5.29b) 



i.e. the dependence slowly approaches the linear func- 
tion. 

To complete this subsection, we consider the effect of 
finite temperature. As for the usual BCS mean field, 
Eq. (5.23) gives the temperature dependence of the width 
of the mean-field gap Ao(T) which contains the scale 
Ao(0) only: 



In 



Tmf 
T 



V- 1 1 

y r =0, 



1 



+ 2) + 477^ 



(5.30) 

where Tc is related to Ao(0) by the usual weak coupling 
BCS relation 



MF 



ne 



^Ao(O) « 1.76Ao(0). 



(5.31) 



and C is the Euler constant. 

Equations (5.31) and (5.30) enable us to estimate the 
upper bound for the mean field crossover temperature. 
For = 4 we find 



-MF 



< IQ-^B, 



(5.32) 



where B is the Zeeman splitting. As the electron g-factor 
in graphene equals to 2, we estimate for the parallel mag- 
netic field B ~ 40T, Tmf — QOmK, which does not seem 
to be non-realistic. The other experimental realization 
could be putting the appropriate insulating ferromagnet 
on the top of the graphene film, so that the Zeeman split- 
ting is caused by the corresponding exchange fields. The 
effective Zeeman splitting in this case may reach thou- 
sands of Kelvin-s^^. 

In the vicinity of the mean-field crossover Tc — T <C Tc, 
we obtain from Eq. (5.30) 



7C(3) 

and C(a;) is the Riemann ^-function. 



Tc{T - T,) « 9.38 Tmf{T - Tmf), 

(5.33) 



D. Microscopic calculation of the coefficients in the 
Free energy and effective action. 

Using the symmetry arguments of Sec. Ill or by the ex- 
plicit calculation, one finds that the solution of the form 
(5.25) is not unique, and, in fact any order parameter of 
the form 



A(n) - ® Q {Ao{T)V{n) + A,{T) 



1 



Vin)] } 
(5.34) 

solves Eq. (5.23), i.e. the fiuctuation effects, studied in 
Sec. IV, are important. Here 4x4 matrix Q is given by 
Eqs. (3.2). 

The free energy and the effective action describes the 
contribution of configurations of the order parameter 
slowly varying in time and space. The time and space 
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gradient terms describes the cost of creating an inhomo- 
geneous configuration, whereas the anisotropy leads the 
certain modes to become massive. To calculate the lat- 
ter ones it is sufficient to consider the homogeneous and 
time independent configurations of the order parameter 
Q, see Eqs. (3.2), whereas to find the former, one has to 
expand in small gradients of Q.^ 



1. Stiffness. 

To find the gradient terms we expand the order pa- 
rameter as 



Q(r) = Q 

6Q{v 



1-\5q\v) 



5Q{v)- 



_fk_ 
(2^ 



(5.35) 



In order to preserve the constraints (3.2b) up to the 
second order perturbation theory in SQ we require 



Q-Qt. QSQ + SQQ^O; 



-SQ. 



(5.36) 



5o = 




FIG. 16: Microscopic calculation of the stiffness and the col- 
lective modes velocity v, in the effective action. Notation is 
defined in Eqs. (5.20) and (5.21). 

Calculating the diagram on Fig. 16, we obtain for the 



^ In principle one can perform the gauge transformation Q — ► 
UQU^ = Sx, and expand the fermionic determinant in terms of 
the non-abehan scalar, U^OtU, and vector, (7^9^ [/, potentials. 
We have chosen not to do it here to avoid careful considerations 
of the terms arising from the non-gauge invariant integration 
cut-off in Eq. (5.22). 



isotropic part of the Free energy 



(2^)2 



Tr^Q(e„i, A: + q)SQqQQ{em,k)SQ _q 
-TvQQ{eyn,k)5QqQQ{ejn,k)5Q_q , 



(5.37) 



where are the fermionic Matsubara frequencies. The 
term in the last line is local and it comes from the in- 
teraction part of the Hamiltonian. Its explicit calcula- 
tion is not necessary because (5(5w=o,fc=o corresponds to 
the motion along the degenerate manifold which can not 
produce any contribution to the Free energy - this re- 
quirement fixes the term unambiguously. 

Expanding the Green functions in powers of small mo- 
mentum q we find 

1 f cpQ f <Pk 



"9 rpY^ 

(2^ e 



(2^)2 



[g2 ^ TidkM^^n,k)5Q,dkM^^,k)5Q- 



Then, the simple power counting shows that the inte- 
gral is contributed by e,vk <C B and therefore the in- 
tegration rule (5.22) may be used. Finally, using the 
explicit expression for the Green function from Fig. 16 
and Eq. (5.36), we obtain after simple algebra 



Fo = 



Pk{T) 



Tr / (fr 



VQ 



Pk{T) 



-f 



(5.38) 



where 



(5.39) 



and C,{x) is the Riemann ('-function. 

As we have already explained, keeping the imaginary 
time derivative terms is valid only for A(T) S> T. At 
such low temperatures we find 

B 



PK = 



47r 



(5.40) 



independently of any logarithmic renormalization from 
higher energies. 

In the opposite case, A(T) <C T only constant in time 
fluctuations are important and we obtain So ~ Fo/T, 
where Fo is given by Eq. (4.4a) and the stiffnesses on 
the mean-field correlation length, S^mf see Eq. (4.3), are 
given by 

Pk{£.mf) = Psi^MF) = f^^) . (5.41) 



27r 



Tmf 



where we used Eq. (5.33). 
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2. Leading anisotropies 



find 



The quadratic anisotropics (4.4b) arises both due to 
the warping of the single electron spectrum (2.13) and 
the short-range interactions (2.14). The diagrammatic 
calculation of these anisotropies is straightforward and it 
is shown on Fig. 17 a). 







+ iny) + iEj, + h.c. 



a) Fx = 



b) Wa = 




FIG. 17: Microscopic calculation of the anisotropic terms 
in the free energy. Diagrams containing largest power of 
ln(Z3/Ao) are marked by star. 



We find F_L = 



r(2) 



where the first contribution 



is due to the warping 



-TrQA.QA^ 



fp 



16t:v^{Rb) \2nT 



Ao 



(5.42a) 



where fp is given by Eq. (5.39). 

The second diagram of Fig. 17 b) describes the effect of 
the short range interaction (2.14). We take into account 
only the most relevant terms as specified in Eq. (5.12b). 
Term proportional to do not contain the inter-valley 
A-matriccs and, thus, do not cause the anisotropy. We 



X Rbv{Rb) 



F1{Rb)Ti- ( QA,y+ F^IRb) ^Tr [QA, 

tJ.=x,y 

Bd^ Ao(T) 



27rv^RB)el+e + ^o(Tr 

(5.42b) 



The logarithmic integral in Eq. (5.42b) is eliminated us- 
ing Eq. (5.26), and we obtain Eq. (4.4b) at the scale 
''0 = £,MF, see Eq. (4.3), with the couplings 



(5.43a) 




(5.43b) 



where fA{x) is defined in Eq. (5.28). 

Equations (5.43) give complete expressions for the 
anisotropies in the free energy (4.4b) in terms of the mi- 
croscopic coupling constants defined on the scale Rb- 
Equations (5.9), (5.12), (5.13), (5.4), (5.18), (5.19), 
(5.28), and (5.29) enable us to find the magnetic field 
dependence of the anisotropy constants. In the case of 
the strong magnetic field B > Bq, see Eq. (5.17), it yields 



-Bixim 



68 B 
12.9Ff (7^)| 



3.9 



Alio) ^ /Ao(T) 



A^(r) 



2ttT 



rii_ w ~12.9BF^{n). (5.44a) 
For the weak magnetic field. B < Bq, we find 




exp 



' TT [In |i + 6.28] 
2 In [In + 6.28] 



0.15 In 



-13/4 



12.9Ff(7^)( 



-12.9BF^{TZ) 



0.151n(|^ 



0.15 In 



(5.44b) 



Comparison of the anisotropies from Eqs. (5.44) with 
the expressions for stiffness (5.41) together with the es- 
timates _Ff ' (7?.) ~ 1 indicates that the weak anisotropy 
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case is possible only for B ^ Bq and at the strong mag- 
netic field the anisotropy is dominating already at the 
mean- field correlation length, see Sec. IV C. 

The very peculiar situation arises if F^{TZ) < 0, as 
we explained before, the value and the sign of this con- 
stant is determined by the details at the distances of 
the order of the lattice constant. In this case the ex- 
change and the warping produce the contribution of the 
different sign and intersection of the Heisenberg line, 
r]z > 0, \r]±\ = see Sec. (IVB2) becomes possi- 

ble. Using Eq. (5.44a), and estimates (5.10) and (5.13), 
we find the estimate for such field Bh 



Bh ^ 10-^Bo 



\Fiin) 



|i^^(7^)Hy 



l^lOSo > Bq. 



XI (7^) 

(5.45) 

The manifestation of this line on the phase diagram will 
be considered in Sec. VE. 

Calculation of the other anisotropics is self-explanatory 
from diagrams Fig. 17 b-c). Though, formally, they are 
of the same order in perturbation theory, they still can 
be classified in powers of ln(S/A) ^ 1. Taking into 
account only the leading logarithmic contributions, and 
reexpressing the logarithmic expression using the mean- 
field equation (5.26), we find for anisotropy coefficients 
in Eqs. (4.4c) - (4.4d) 



, B N[F-{n)-iJim\ \ /Ao(T) 



(5.46a) 



9 js 

a^MF)^ — [T-{n)fD 

TT 



N 



' Ao(0) ' 



^/a (life)) 

Ao(r) - 

2ttT 



(5.46b) 



where D{x) = xe ^ , the functions fp^A are defined 
in Eqs. (5.39) and (5.28). Deriving Eq. (5.46a), we 
took into account that k is important only if 77^ — > 0, 
see Sec. IV Bl. Thus, the irrelevant constants, see 
Eqs. (5.12a) had to be taken into account. 

Expressions (5.46) together with Eqs. (5.12a) and 
(5.15) shows that those anisotropics are very small, so 
we will not analyse their asymptotics further. 



E. Phase diagram in B-T plane. 

This subsection combines the symmetry analysis of 
Sec. IV with the microscopic calculation of the Free en- 
ergy couplings in Sec. VD. As the result, we will con- 
struct the phase diagram in the plane determined by the 
Zeeman splitting B and by the system temperature, T. 
As we have discussed in Sec. IV, the interesting phase 
transitions are determined by the thermal fluctuations 



and it is convenient to introduce dimensionless Ginzburg- 
Levanyuk parameter 



Gi 



B 



«1, 



(5.47) 



characterising strength of such fluctuations. Here the 
mean-field transition temperature is given by Eqs. (5.31) 
and (5.27) - (5.29). 

Apparently, not all the regions of the phase diagram 
of Fig. 9 can be explored by varying B, T and we will 
restrict ourselves with two most realistic, as we believe, 
cases. Namely, we assume that the absolute values of 
the interaction constants |Ff'^(a)| ~ 1. For the sake of 
concreteness, we assume J^i"(a) > 0. [Case of F^{a) < 
is obtained by the replacement A''i 7V4.] 

Let us consider first the case of Ff < 0. Then, ac- 
cording to Eqs. (5.43), one finds ry^.x < 0. Mean-field 
diagram obtained from Fig. 4 is trivial and includes the 
continuous transition from the disordered normal state, 
to the spin- flux state. Fig. 3 (d), of the excitonic insula- 
tor, see Fig. 18. 

The fine structure of the phase diagram. Fig. 18 (b), 
is obtained from the general Fig. 9, by using the micro- 
scopic expression for the Free energy couplings derived 
in Sec. VD. The topological structure of the phase dia- 
gram is most easily obtained by the mapping of the paths 
in B,T plane to the path in {tIztVI-) plane as shown in 
Fig. 18 (c). 

The positions of the transitions lines on the phase di- 
agrams are obtained by combining the phenomenolog- 
ical results of Sec. IV and the microscopic analysis of 
Sec. VD. 

For instance, using Eqs. (4.15) and (5.41), we obtain 
in the limit of the small vortex fugacity 



MF 



T 



TmfGi 



hne (i); 
line (ii). 



(5.48a) 



Analogously, using Eqs. (4.20), (4.14) and (5.44b) we ob- 
tain a position of the Ising line. With the logarithmic 
accuracy, we find for line (iii) of Fig. 18 (b): 



MF 



T 



TmfGi 



iln^; Bocxp{-l/Gi) <B <s:BoGi 
^0; B> BoGi. 



(5.48b) 

Case of Ff < is more sophisticated. As we noticed 
in Sec. VD, coefficient r/^ changes its sign as the func- 
tion of the magnetic field and at some point crosses the 
Heisenberg line at field Bh- At the mean-field level, it 
corresponds to the continuous transition between two- 
kinds of excitonic insulator: spin flux state, see Fig. 3 
(d). and the link centered spin density wave, see Fig. 3 
(b,c). ^ 

Similarly to the previous case, the fine structure of the 
phase diagram. Fig. 19 (b), is obtained from the general 
Fig. 9, by the mapping of the paths in B, T plane to the 
path in (77^, ri±^ plane as shown in Fig. 19 (c). 
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FIG. 18: Phase diagram of the graphene in the parallel mag- 
netic field for the short range interaction constant Fl > 0. 
a) Mean field structure of the phase diagram; b) The "fine" 
structure of the phase diagram in the close vicinity of the 
mean-field transition temperature, Tmf{B); c) Relation of 
the phase diagram (b) to the more phenomenological phase 
diagram of Fig. 9. 



FIG. 19: Phase diagram of the graphene in the parallel mag- 
netic field for the short range interaction constant Fz < 0. 
a) Mean field structure of the phase diagram; b) The "fine" 
structure of the phase diagram in the close vicinity of the 
mean-field transition temperature, Tmf{B); c) Relation of 
the phase diagram (b) to the more general phenomenological 
phase diagram of Fig. 9. Cross-section denoted by a, 13, 7 are 
shown on Fig. 18 c). 



The Berezinskii-Kostcrlitz-Thouless transition lines 
(i), (ii) on Fig. 19 (b) are still determined by the ex- 
pressions (5.48a). Ising line (iv) Berezinskii-Kosterlitz- 
Thouless line and the are found from Eqs. (4.40) and 
(5.41). Expanding 



A?; ; 



B-B 



H 



Bh 



we obtain with the logarithmic accuracy for 
Bh exp(-G^) < |6 - 6^1 < Bh: 



1] r (■ ^ 

— til— — — ; line (iv); 



Tmf^T I 8 B-Bh 
TmfGi I I 



(5. 
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For larger B the Ising line (iv) approaches the mean- 
field temperature. The Berezinskii-Kosterlitz-Thouless 
transitions lines (v), for B <C Bh can be found us- 
ing Eqs. (4.20), (4.14), and (5.44b). For the fields 
Bo exp {~l/Gi) < S < BoGi this yields 

~^ 7^ ~ Ime v ; 5.49b 

ImfGi 4 B 

The Berezinskii-Kosterlitz-Thouless line (vi) turns out to 
lie outside the fiuctuation region due to the large numer- 
ical factor in Eq. (4.44) and smallness of C, in Eq. (5.46b). 
We will not write-down its asymptotic behaviour. 

This completes our analysis of the structure of the 
phase diagram of graphene in the parallel magnetic field, 
characterised by the Zeeman splitting B. 

VI. SUMMARY AND CONCLUSIONS 

In this paper we have discussed two problems concern- 
ing clean graphene: (i) possible effects of in-plane mag- 
netic field in facilitating a formation of excitonic conden- 
sate and (ii) a role of the long range Coulomb interac- 
tion and its influence on other interactions in the system. 
The second topic is more general than the first and has 
a broader significance, though, as far the paper goes, it 
was discussed in the second part. 

In zero magnetic field graphene is a gapless semicon- 
ductor with two Fermi points in the Brillouin zone (val- 
leys). In-plane magnetic field pushes up- and down-spin 
bands in opposite directions transforming the system into 
a metal with extended Fermi surfaces for electrons and 
holes of opposite spin. There are two such Fermi sur- 
faces corresponding to two possible valley indices. Elec- 
trons and holes attract through the Coulomb interaction 
which creates a possibility of exciton condensation along 
the lines first described by Keldysh and Kopaev^^. The 
maximal possible symmetry of the order parameter is 
U(2), lattice effects bring it down to U(l). The system 
in its low temperature phase is an insulator with a gapless 
collective mode corresponding to fiuctuations of spin den- 
sity in the directions transverse to the applied magnetic 
field. The above is a brief summary of the discussion of 
Sections II-IV. Section IV also contains a detailed phase 
diagram. The discussion in these sections dealt with the 
Landau-Ginzburg free energy functional written purely 
on symmetry grounds where various energy scales of the 
system enter as parameters. 

Section V contains a microscopic analysis tailored es- 
pecially for graphene. The ultimate goal was to obtain 
estimates for the critical temperature and various pa- 
rameters of the phase diagram, but a byproduct of the 
analysis is a study of how the strongest interaction in 
graphene - the Coulomb interaction affects the spectrum 
and renormalizes other interactions (such as the short 
range exchange). As is well known, the Coulomb inter- 
action in graphene, measured by its dimensionless value 
g = ire^ /2v is quite strong at energies of the order of 



the bandwidth. We have found, however, that the ef- 
fective coupling steadily diminishes at low energies and 
asymptotically vanishes at £" = 0. The scale dependence 
of g is rather slow and is given by Eq.(5.4). This for- 
mula includes an important scale TZ, which we estimate 
for graphene as being of order of 10^ — 10'^ lattice con- 
stants. This scale separates the region of relatively strong 
interaction where g{r) decreases as a power law, from the 
region of weak coupling where g{r) ^ [In r] ~^ . It also sets 
the scale Bo for the magnetic field (5.17)[our estimate is 
6o « 10 - IGGK]. 

The renormalization process is drastically altered at 
energies of order of the applied magnetic field B. The 
field sets the ultraviolet cut-off for the physics of excitonic 
insulator. However, the upper cut-off for its collective 
excitations is much lower and is set by the value of the 
mean field gap Aq. One may anticipate that the latter 
energy scale is exponentially small in comparison with 
the cut-off B. This is indeed the case, but fortunately 
the inverse coupling constant l/g{r ^ B^^), which stays 
in the exponent, depends rather weakly on the magnetic 
field [see (5. 29a, 5. 29b)] so that the magnitude of Aq is not 
that small. Our estimate is that in fields of the order of or 
stronger than lOT the mean field temperature is Ti\ip < 
lO^^B. This makes it possible to observe the excitonic 
effects described in this paper in the temperature range 
of tens of mK (see more discussion in Section V C). 

Though the long range Coulomb interaction is cer- 
tainly the main player, its renormalization drags with 
itself weaker interactions, such as the short range ex- 
change, strengthening them at low energies. Such inter- 
actions break the U(2) symmetry present at low energies 
when only the long range Coulomb interaction is taken 
into account. The analysis of Section V C.2 demonstrates 
that the U(2) symmetry survives in the excitonic insu- 
lator only at B << Bo when the estimated transition 
temperatures are probably too low for the effect to be 
observed. In the realistic region B > Bq the anisotropy is 
strong. The expected phase diagrams in B — T plane are 
given on Figs. 18,19 [they differ by a sign of a certain in- 
teraction parameter which on the current stage remains 
unknown]. Strong fields of order of lOT or more, which 
are required to make the excitonic insulator observable 
at temperatures above tens of mK, will probably put one 
in the regime marked by /3 or 7 on the phase diagram 
Figs. 18, 19. The spin configuration corresponding to 
this regime is depicted on Fig. 3 d) and corresponds to 
the spin-fiux phase. 
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FIG. 20: Cuts on the x — y plane attached to each topological 
defect in the five- vortex configuration. 



must be multi-valued function of the coordinate. 

In order to avoid the consideration of the multi- valued 
field we introduce cuts parallel to, say, x axis connecting 
each vortex or half- vortex with the boundary of the sys- 
tem [any physical quantity, obviously, does not depend 
on the choice of the cut], see Fig. 20, and consider all 
the matrices to be single valued function but the phases 
experiencing the discontinuity on the cuts. 

The current (A2) should be continuous thus we modify 
the definition as 



APPENDIX A: DERIVATION OF EQ. (4.11). 

Let us re-write Eq. (4.4a') in a form 



where 



TP PK ^ -.2 , Ps - Pk 



(Al) 



(A2) 



The topological defects (4.7) and (4.9) are determined by 
the condition 



I.e. 



-Tr In V 



(A3) 



(A4) 



A/-(l/2) 



(A5) 



where A/"'^' and M^^l"^^ are the number of the vortices 
and the vertices respectively, ccj, j/j are their coordinates, 
Ij = ±1 is the corresponding vorticity and nij is the unit 
vector characterising the spin of the jth half-vortex. 

Writing the summation over the vortex and half- vortex 
coordinates explicitly, wc obtain 



E 



A/-(l/2) ^(1/2) 

n j: 



^"(1/2) I 

A/'(i/2)=o ■ J = l ;(.i/2)=±i 



dxf/'^dyf/'Um, 



VV 



X exp 



dxdy 



4T 8T ^ 



(A6) 



where the matrix current is given by Eq. (A5). 

After introducing the dual 2x2 matrix field h ^ h\ Eq. (A6) acquires the form 



X exp 



dxdy 



dx^j 'dy)j 



E 



^(1/2) ^(1/2) 

n E 



/V(l/2)l 

Ar(l/2)=0 ■ j=l ;(l/2)^±x 



dxf'^Uyf'^'^drr., 



VVVh 
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4T ' 4pK 
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8T 



8ps 8pK 



Trdxhj +iTrjydxh 



(A7) 
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Substituting Eq. (A5) into Eq. (A7); integrating the terms with S functions by parts and summing over Ij, we find 



VVVh E /^^2cos.T.%W) 



j=i 

^(1/2) j^(l/2) (1/2) , (l/2)dm. 



X 



^"(1/2) I 11 



X exp 



47rr§ 



-2 COSTT 



Tr/i(' 



,(1/2). 1 + 



PA- 



ST 



T T \ 
8ps 8pK J 



Trd^h) +TrV'<dyVda:h 



(A8) 



After the summation over A/"^^/^-* , A^^'^^ , and integration Then Eq. (B3) becomes 
over m.j, we obtain Eq. (4.11). 

U = -{ladxL - rad^ra) + T{liri){l2r2) 

APPENDIX B: ANALYSIS OF THE ISING / 

PHASE TRANSITION. ^{v±/T + P^aFiW + i{v±/T - PM)r2h (B6) 



In the vicinity of T = i:pK, where the mutually dual 
cosines have the same scaling dimension 1, partition 
(4.25) can be mapped to the quantum many-body prob- 
lem at zero temperature and then refermionized. The 
reader can consult Ref. 23 where the necessary infor- 
mation about 2D Ising model is provided. Choosing y 
coordinate for imaginary time, we re-write the classical 
Eq. (4.25) 



Z oc Tr exp —Ly 



n 



dxH 
T 



+ cos 20 - 



2pK 



^2 ^"^M 

i6{x — x') 



(Bl) 



dj{x)-e{x') 

Then using the fermionization rules 
R(x) (X exp 
L{x) cx exp —i(j){x) + iirO 



(B2) 



we write down the corresponding ID quantum fermionic 
Hamiltonian density as 

n = i{L+dxL - R+dxR) + [t^T/pk - l)R+RL+L 
+T]^{R+ L + L+ R) + p3i{R+ L+ + LR) (B3) 

It is convenient to decompose the Dirac spinor into the 
real (Majorana) components: 

R = ri+ir2, L = k + ih (B4) 
where the corresponding operators are real (r^ = r+ , = 
Z+) and satisfy the following commutation relations: 

{raixi),rb{x2)} = 5abS{xi2), (B5) 
{la{xi)Jbix2)} = SabS{xi2),{ra{xi),lb{x2)} = 



where r = [T/{ttpk) — 1]. Hamiltonian (B6) describes 
two quantum Ising models coupled by the energy density 
operators. The original order parameter operator can be 
expressed in terms of order and disorder parameters of 
the Ising models a and p, see Ref. 23 for the correspond- 
ing definitions: 



e"^ = (Ti(T2 + ipiP2 



(B7) 



The sign of Majorana mass in the Ising model plays 
an important role determining what operator (cr or p) 
acquires a vacuum expectation value. Then from (B7) it 
is clear that this operator acquires a finite expectation 
value when the masses of the two species of Majorana 
fermions have the same sign so that either cr^ or pa fields 
condense simultaneously. The high temperature phase is 
characterised by masses of different sign. 

At small |r| << 1 we can use perturbation theory to 
write the equations for the masses (let us choose ri± > 0): 

mi = ii]±/T + P34) > 0; 

m2 = {V±/T - P34) + ln(l/i?*|mi|); (B8) 

It follows that 7712 changes sign at 



77_L /T — /i34 

Tc/npK - 1 = — In 



'n±/T + P34 \R*i'n± + P34) 



(B9) 



where the second order phase transition from disordered 
(high temperature) to the ordered (low temperature) 
state takes place. Equation (B9) agrees with Eq. (4.29). 
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